/*
 * ibdscanner.cpp
 *
 *  Created on: Sep 27, 2011
 *      Author: lindenb
 */
#include <iostream>
#include <algorithm>
#include <cerrno>
#include <cmath>
#include <cassert>
#include <fstream>
#include <cstdlib>
#include <set>
#include <vector>
#include <cstdlib>
#include <cstring>
#include <cmath>
#include <limits>
#include "xstdio.h"
#include "xstring.h"
#include "xzlib.h"
#include "throw.h"
using namespace std;

#define WHERE(a) cerr << __FILE__ <<":"<< __LINE__<<":" <<  a << endl
#define DEBUG_LINE(s) s->trace_line=__LINE__
class Segment
	    {
	    public:
			std::set<string*> names;
			std::string chrom;
			long start;
			long end;
			double p_value;
			int trace_line;
			void print(std::ostream& out) const
				{
				out
						// << "#"<< trace_line << "\t"
						<< chrom << "\t" << start << "\t"<< end << "\t"<< p_value << "\t"
						<< names.size() << "\t"
						;
				bool first=true;
				for(std::set<string*>::const_iterator r=names.begin();r!=names.end();++r)
					{
					if(!first) out << ";";
					first=false;
					out << *(*r);
					}
				}
	    };

std::ostream& operator<<(std::ostream& out,const Segment& cp)
	{
	cp.print(out);
	return out;
	}

class Scanner
    {
    public:
		set<std::string> names;
		ostream* out;
		double max_diff;


	static bool compareseg(const Segment* a,const Segment* b)
		{
		if(a->names.size()> b->names.size()) return false;
		if(a->names.size()< b->names.size()) return true;
		return (a->end-a->start)> (b->end-b->start);
		}



	vector<Segment*> segments;


	Scanner():out(&cout),max_diff(0.1)
	    {
	    }

	~Scanner()
		{

		}


#define INSERT_NAMES_TO(A,B) A->names.insert(B->names.begin(),B->names.end())


	bool insert(Segment* seg,int depth)
		{
		bool inserted=false;
		vector<Segment*> buffer;

		//if(seg->end==159288) WHERE(depth << " "<< segments.size() );
		assert(seg->start<=seg->end);
		assert(seg->start>=0);

		for(size_t i=0;i< segments.size();++i)
			{
			Segment* other=segments[i];
			assert(other->start<=other->end);
			if(other->chrom.compare(seg->chrom)!=0 ||
				fabs(other->p_value - seg->p_value)> max_diff ||
				other->start>seg->end ||
				other->end< seg->start
				)
				{
				continue;
				}

			inserted=true;
			if(other->start==seg->start && other->end==seg->end)
				{
				INSERT_NAMES_TO(other,seg);
				break;
				}
			else if(other->start<= seg->start &&
					other->end >= seg->end)
				{

				//left
				if(other->start< seg->start)
					{
					Segment* newseg=new Segment;
					DEBUG_LINE(newseg);
					newseg->chrom=seg->chrom;
					INSERT_NAMES_TO(newseg,seg);
					newseg->start=other->start;
					newseg->end=seg->start - 1;
					newseg->p_value=other->p_value;

					buffer.push_back(newseg);
					}
				//right
				if(other->end> seg->end)
					{
					Segment* newseg=new Segment;
					DEBUG_LINE(newseg);
					newseg->chrom=seg->chrom;
					INSERT_NAMES_TO(newseg,seg);
					newseg->start = seg->end +1;
					newseg->end = other->end;
					newseg->p_value=other->p_value;

					buffer.push_back(newseg);
					}
				//mid
				other->start=seg->start;
				other->end=seg->end;

				INSERT_NAMES_TO(other,seg);
				}
			else if(other->start >= seg->start && other->end <= seg->end)
				{
				if(other->start> seg->start)
					{
					Segment* newseg=new Segment;
					DEBUG_LINE(newseg);
					newseg->chrom=seg->chrom;
					INSERT_NAMES_TO(newseg,seg);
					newseg->start = seg->start;
					newseg->end = other->start -1;
					newseg->p_value= seg->p_value;

					buffer.push_back(newseg);
					}
				if(other->end< seg->end)
					{
					Segment* newseg=new Segment;
					DEBUG_LINE(newseg);
					newseg->chrom=seg->chrom;
					INSERT_NAMES_TO(newseg,seg);
					newseg->start = other->end +1;
					newseg->end = seg->end;
					newseg->p_value= seg->p_value;
					buffer.push_back(newseg);
					}
				//mid
				//other->start=other->start;
				//other->end=other->end;
				INSERT_NAMES_TO(other,seg);
				}
			else if(other->start>= seg->start && other->end >= seg->end)
				{
				if(other->start> seg->start)
					{
					Segment* newseg=new Segment;
					DEBUG_LINE(newseg);
					newseg->chrom=seg->chrom;
					INSERT_NAMES_TO(newseg,seg);
					newseg->start = seg->start;
					newseg->end = other->start -1;
					newseg->p_value= seg->p_value;
					buffer.push_back(newseg);
					}
				if(other->end> seg->end)
					{
					Segment* newseg=new Segment;
					DEBUG_LINE(newseg);
					newseg->chrom=seg->chrom;
					INSERT_NAMES_TO(newseg,seg);
					newseg->start = seg->end +1;
					newseg->end = other->end;
					newseg->p_value=other->p_value;
					buffer.push_back(newseg);
					}
				INSERT_NAMES_TO(other,seg);
				other->end=seg->end;
				}
			else if(other->start<= seg->start && other->end <= seg->end)
				{
				if(other->start<seg->start)
					{
					Segment* newseg=new Segment;
					DEBUG_LINE(newseg);
					newseg->chrom=seg->chrom;
					INSERT_NAMES_TO(newseg,seg);
					newseg->start=other->start;
					newseg->end=seg->start -1;
					newseg->p_value=other->p_value;
					buffer.push_back(newseg);
					}
				if(other->end< seg->end)
					{
					Segment* newseg=new Segment;
					DEBUG_LINE(newseg);
					newseg->chrom=seg->chrom;
					INSERT_NAMES_TO(newseg,seg);
					newseg->start = other->end +1;
					newseg->end = seg->end;
					newseg->p_value= seg->p_value;
					buffer.push_back(newseg);
					}
				INSERT_NAMES_TO(other,seg);
				other->start=seg->start;
				}
			else
				{
				*out << "other:";
				other->print(*out);
				*out << "seg:";
				seg->print(*out);
				THROW("boum");
				}
			}


		if(!inserted)
			{
			//WHERE(segments.size() << " "<< depth <<" " << *seg);
			segments.push_back(seg);
			}
		for(size_t i=0;i< buffer.size();++i)
			{
			if(buffer[i]->end-buffer[i]->start > 1 )
				{
				//WHERE(*buffer[i]);
				insert(buffer[i],depth+1);
				}
			}
		return true;
		}

	void run(Reader* reader)
	    {
	    vector<string> tokens;
	    vector<string> headers;
	    char* line;

	    if((line=reader->readLine())==NULL) THROW("Cannot read first line");
	    if(line[0]!='#') THROW("Expected first line to start with #");
	    XCString::split(line,"\t",headers);

	    while((line=reader->readLine())!=NULL)
			{
			if(line[0]==0 || line[0]=='#') continue;

			XCString::split(line,"\t",tokens);
			if(tokens.size()!=headers.size()) THROW("Runtime error number of column!=header in "<< line);
			pair<set<string>::iterator,bool> p=names.insert(tokens[0]);
			//string& pair_name2=*(p.first);
			const string* pair_name=&(*(p.first));
			char* p2=NULL;
			Segment* seg=new Segment;
			DEBUG_LINE(seg);
			seg->names.insert((string*)pair_name);
			seg->chrom=tokens[1];
			seg->start=atoi(tokens[2].c_str());
			seg->end=atoi(tokens[3].c_str());
			seg->p_value=strtod(tokens[5].c_str(),&p2);
			insert(seg,0);
			}
	    sort(segments.begin(),segments.end(),Scanner::compareseg);
		for(size_t i=0;i< segments.size();++i)
			{
			//if(segments[i]->names.size()<2) continue;
			*out << *(segments[i]) << endl;
			}

	    }
    };

static void usage(std::ostream& out,int argc,char** argv)
		{
		out << argv[0]<< "\nLast compilation:"<< __FILE__ << " "<< __LINE__
		<<"\nPierre Lindenbaum PhD.\n"
		<<"Options:.\n"
		<<" -o (fileout)\n"
		;
		out << endl;
		}

int main(int argc,char** argv)
    {
    Scanner scanner;
  
    int optind=1;

    char* fileout=NULL;
    fstream* fout=NULL;
    while(optind < argc)
	    {
	    if(std::strcmp(argv[optind],"-h")==0)
			{
			usage(cout,argc,argv);
			return EXIT_SUCCESS;
			}
	    else if(std::strcmp(argv[optind],"-o")==0 && optind+1<argc)
			{
			fileout=argv[++optind];
			}
	    else if(std::strcmp(argv[optind],"--diff")==0 && optind+1<argc)
			{
			scanner.max_diff=atof(argv[optind]);
			}
	    else if(std::strcmp(argv[optind],"--")==0)
			{
			++optind;
			break;
			}
	    else if(argv[optind][0]=='-')
			{
			std::cerr << "unknown option " << argv[optind] << std::endl;
			usage(cout,argc,argv);
			return (EXIT_FAILURE);
			}
	    else
			{
			break;
			}
	    ++optind;
	    }


    if(fileout!=NULL)
	    {
	   fout=new fstream(fileout,ios::out);
	    if(!fout->is_open())
		    {
		    cerr << "Cannot open \""<< fileout << "\" " << strerror(errno) << endl;
		    return EXIT_FAILURE;
		    }
	    scanner.out=fout;
	    }





    try
	    {
	    if(optind==argc)
			{
			GZReader input;
			scanner.run(&input);
			}
	    else if(optind+1==argc)
			{
			GZReader input(argv[optind++]);
			scanner.run(&input);
			input.close();
			}
	    else
			{
			cerr << "Illegal number of arguments." << endl;
			return EXIT_FAILURE;
			}
	    }
    catch(std::exception& err)
	    {
	    cerr << "An error occured:" << err.what()<<  endl;
	    }
    catch(...)
	    {
	    cerr << "An error occured" << endl;
	    }
    if(fout!=NULL)
	    {
	    fout->flush();
	    fout->close();
	    scanner.out=&cout;
	    }
    return EXIT_SUCCESS;
    }
